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ABSTRACT 

A  two-dimensional  flow  governed  by  the  incompressible  Navier- Stokes  equations  with  a 
steady  spatially  periodic  forcing  (known  as  the  Kolmogorov  flow)  is  numerically  simulated. 
The  behavior  of  the  flow  and  its  transition  states  as  the  Reynolds  number  Re  varies  is  inves¬ 
tigated  in  detail,  as  well  as  a  number  of  the  flow  features.  A  sequence  of  bifurcations  is  shown 
to  take  place  in  the  flow  as  Re  varied.  Two  main  regimes  of  the  flow  have  been  observed: 
small  and  large  scale  structure  regimes  corresponding  to  different  ranges  of  Re.  Each  of 
the  regimes  includes  a  number  of  quasiperiodic,  chaotic  and  relaminarization  windows.  In 
addition,  each  range  contains  a  chaotic  window  with  non-ergodic  chaotic  attractors.  Spa¬ 
tially  disordered,  but  temporally  steady  states  have  been  discovered  in  large  scale  structure 
regime.  Features  of  the  diverse  cases  are  displayed  in  terms  of  the  temporal  power  spectrum, 
Poincare  sections  and,  where  possible,  Lyapunov  exponents  and  Kaplan- Yorke  dimension. 
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1  Introduction 


It  is  now  well  established  that  a  number  of  dissipative  chaotic  fluid  flows 
may  require  relatively  few  dimensions  for  their  description1-4.  This  has  boon 
examined  perhaps  in  greatest  detail  for  fluid  models  such  as  the  Ginzburg- 
Landau  (G-L)  equation5-8,  an  equation  which  results  from  the  study  of  the 
critical  point  in  a  variety  of  stability  problems9-12.  For  closed  fluid  systems 
such  as  Raylcigh-Bcnard  (R-B)  convection  and  Taylor-Couette  flow  there  is 
abundant  evidence  that  the  early  stages  of  chaos  are  low  dimensional.  For 
R-B  convection  there  also  exists  a  numerical  simulation  which  details  the 
variation  of  dimension  with  Rayleigh  number18. 

In  this  paper  we  consider  a  simple  fluid  flow  which  exhibits  chaotic  behav¬ 
ior  and  which  is  amenable  to  detailed  examination.  'Phis  flow,  governed  by 
the  incompressible  Navier-Stokes  equations,  is  generated  and  maintained  by 
a  spatially  periodic  time  independent  forcing  of  the  fluid.  It  is  sometimes  re¬ 
ferred  to  as  the  Kolmogorov  flow  (Arnold  and  Meshalkin14).  The  stability  of 
the  primary  flow  has  been  reported  on  by  Meshalkin  and  Sinai15  and  Green11'. 
That  flow  is  known  to  become  unstable  at  a  crit  ical  Reynolds  number.  /iY,.,  of 
v2  beyond  which  a  st  at  ionary  cellular  pat  t eru  appears1 ' .  At  suliicient  ly 

high  !\<  the  cellular  pattern  itself  becomes  unstable  and  the  resulting  flow 
is  unsteady  and  chaotic.  Kraiehnan18  and  Batchelor111  have  shown  that  two 


dimensional  turbulent  flows  exhibit  an  energy  cascade  toward  low  wavenum¬ 
bers  (the  infamous  reverse  cascade)  and  an  enstrophy  cascade  toward  high 
wavenumbers.  This  phenomenon  is  apparent  in  our  calculations.  An  alter¬ 
nate  view  of  this  effect,  in  the  present  instance,  is  that  the  flow  is  unstable 
to  long  wavelength  disturbances.  The  general  instability  of  short  wavelength 
structures  to  long  wavelength  disturbances  has  been  investigated  by  Yakhot 
and  Sivaskinsky,  and  under  sufficiently  anisotropic  conditions,  their  study 
lead  to  the  notion  of  negative  viscosity  for  the  long  wavelengths20,21,  since 
the  growth  rate  is  proportional  to  the  square  of  the  wavenumber. 

Our  investigation  is  in  large  part  numerical  in  nature.  We  examine  in 
some  detail  the  behavior  of  the  flow  and  its  transition  states  for  a  range  of  the 
bifurcation  parameter  Re.  We  use  Poincare  sections  to  examine  the  chaotic 
attractor  of  the  flow,  and  apply  the  Kaplan- Yorke  formula22  to  calculate 
the  dimension  of  the  chaotic  attractor  when  it  is  possible.  A  number  of 
remarkable  features  of  the  flow  were  discovered.  These  include  windows  of 
relaminarization  in  the  parameter  space,  chaotic  windows  with  non-cryodic 
attractors  and  spatially  disordered  but  temporally  steady  states. 

We  will  be  careful  throughout  this  paper  to  use  the  term  chaotic  instead  of 
the  term  turbulent.  In  particular  we  are  referring  to  Hows  which  when  looked 
at  m  the  time  domain  have  broad  band  energy  spectra.  Those  Hows  also 
exhibit  sensitive  dependence  on  initial  data  (they  have  positive  Lyapunov 

•> 


exponents).  Thus  by  generally  accepted  criteria  the  (lows  are  temporally 
chaotic.  However,  as  will  be  seen,  relatively  few  spatial  Fourier  modes  are 
active  and  we  say  we  have  chaos  rather  than  turbulence.  Nevertheless  the 
flow  is  chaotic  in  its  Eulerian  description  and  should  not  be  confused  with 
Lagrangian  turbulence 23-25  in  which  the  velocity  fields  are  to  a  high  degree 
laminar. 


Recently  we  became  aware  of  the  investigation  of  the  Kolmogorov  Flow 
by  She20’27.  While  there  is  some  unavoidable  overlap,  the  two  investigations 
are  largely  complimentary.  She’s  work  is  also  mainly  numerical  and  concen¬ 
trates  on  the  high  end  of  the  Re  range.  (A  minor  difference  is  that  he  uses 
a  higher  value  for  the  frequency  of  the  forcing  function.)  The  mathemati¬ 


cal  tools  used  also  differ.  We  will  make  extensive  use  of  the  power  spectra, 
Poincare  sections,  Lyapunov  spectra  and  flow  visualizations  for  the  investi¬ 
gation,  while  She  relics  mainly  on  time  series  and  phase  portraits.  Several 
issues  are  viewed  differently  by  us  and  we  comment  on  these  in  the  text. 
She’s  further  work  done  in  collaboration  with  Nicolaenko  is  concentrated  on 


2  Kolmogorov  Flow 


We  consider  the  forced  incompressible  Navier-Stokes  equations  in  two  di¬ 
mensions, 


+ 

II 
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II 

o 

(2.1) 

—  u'  +  Vp'/p  =  v  v'2u'  +  X  sin  ky'ex 

(2.2) 

In  preparation  for  the  normalization  given  below,  dimensional  quantities  are 
denoted  by  primes,  er  denotes  the  unit  vector  in  the  or-direction.  This  system 
has  the  solution  or  fixed  point , 

u'  =  U'  =  sin  ky'.  (2.3) 

A  Reynolds  number  for  the  flow  may  be  naturally  based  on  the  maximum 
speed  of  this  field,  \/uk 2  and  the  length  scale  k~ 1 ,  thus 

Re  =  X/(k3’'2).  (2.4) 

(A  Reynolds  could  be  based  on  the  wavelength  and  maximum  speed  in  which 
case  ‘It Re.  is  the  Reynolds  number). 

The  linear  stability  of  (2.3)  can  be  analyzed  by  writing 

u'  =  IT  +  /)u'(y)  exp(m k.r  +  ert)  (2.5) 

and  investigating  the  resulting  Orr-Sommerfeld  equation  gotten  by  lineariza¬ 
tion  in  the  small  pint urbat ion,  It  is  found  that  instability  fust  sets  in 
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Re  =  Rec  =  \/2 


(2.6) 


at15-16 

and  is  due  to  disturbances  of  infinitely  long  wavelength.  Moreover,  it  can  be 
shown  that  disturbances  with  ak  >  1  never  destabilize  U'. 

Thus  according  to  the  stability  analysis,  the  domain  must  be  large  in 
some  sense  in  order  for  us  to  include  growing  modes.  For  computational 
purposes  we  fix  the  domain  size  and  regard  the  spatial  forcing  frequency 
as  the  bifurcation  parameter.  With  this  in  mind  we  choose  the  following 
normalization, 

u  =  uk2 u'/x,  x  =  kx’/n,  t  =  xt' /(nki/),p  =  pi/2k4 /\2,  (2.7) 

where  n  is  an  integer  which  specifies  the  spatial  frequency  of  the  forcing. 
Under  this  normalization  the  Navier-StoKes  equations  become, 

V  ’  u  =  0  (2.8) 

du  1  2  7l2  .  .... 

* +vp=nv  u+nc's"""'  ('-9) 

with 

0  =  nx/(v2k3)  —  v  Rc  (2.10) 

the  bifurcation  parameter.  In  this  format  the  critical  value  of  U  is 

$2,-  =  nRcc  —  n\J 2.  (2.1  1 ) 

r> 


Let  us  now  briefly  discuss  the  symmetry  groups  for  the  solutions  of  these 
equations  which  we  will  refer  to  as  the  Kolmogorov  flows31.  There  are  in  fact 
two  groups,  one  discrete  and  the  other  continuous. 

Referring  back  to  the  notation  of  (2.9),  n  refers  to  the  number  of  vertical 
cycles  in  the  horizontal  forcing.  If  [tz(x,  j/,  t),  v(x,  y,  t)]  =  u  is  a  solution  then 

/u  =  [(-)'■“  ((-)**>!(  +  *4)  ((-)‘x,»  +  40]  k  =  0, .. .  ,2n  -  1, 

(2.12) 

are  also  (independent)  solutions  of  the  problem,  i.e.  they  satisfy  the  equations 
and  conditions  of  the  problem.  The  generator  of  this  cyclic  group,  g ,  is  a 
glide-reflection 32  of  half  wavelength.  In  addition  rotation  through  ~  is  another 
group  generator,  i.e. 

Ru  =  [-u(-x,  -t/),-u(-x,  -y)]  (2.13) 

is  also  a  solution.  In  all  there  are  4 n  elements  in  the  discrete  symmctrv 
group  for  the  Kolmogorov  flow.  Thus  for  the  case  being  considered,  n  —  4, 
we  obtain  a  sixteen-fold  increase  in  the  number  of  possible  solutions  which 
we  can  view. 

The  remaining  symmetry  group  is  the  group  of  translations  in  the  ,r- 
directions,  i.e.  if  u  is  a  solution  then 

//H  --  ["(.r  +  /,'/).  r(.r  4-  /,//)]  (2.14) 

are  also  solutions  to  the  problem  for  all  0  <  /  <  2~. 
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3  Computational  Method 


Based  on  the  normalization  given  in  the  previous  chapter  we  consider  the 
computational  domain 

0<x,y<27r  (3.1) 

with  periodic  boundary  conditions  in  the  two  directions  (so  the  flow  lies  on 
the  torus). 

The  natural  method  of  choice  for  the  investigation  of  (2.8,9)  under  the 
27r-periodic  boundary  conditions  is  the  pseudo-spectral  method  using  Fourier 
expansions33.  All  derivatives  are  evaluated  in  the  spectral  space  and  updates 
are  performed  on  the  expansion  coefficients  of  u(x,  t).  On  the  other  hand  the 
non-linear  terms  are  evaluated  in  the  physical  space  and  then  transformed 
back  to  the  spectral  space.  The  needed  transformations  are  all  efficiently 
done  by  the  means  of  the  Fast  Fourier  Transforms. 

In  a  variation  on  the  usual  formulation  we  have  used  the  vector  expansion 
set 

Biu(x)  =  “ — exp  (mi  •  x),  (3.2) 

where  m  is  the  integer  wavenumber  pair  (?/),.  m2).  Using  such  a  basis  we  can 
expand  the  velocity  field  u  as 

u(x./)  X>m(0Hm(x)  (3.3) 

111 
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where  the  sum  is  performed  over  all  the  components  of  the  wave  vector  m. 

The  basis  elements  are  easily  seen  to  be  divergence  free, 

V-Bm  =  0  (3.4) 

and  complete  in  the  space  of  square  integrable  2~-periodic  divergence- free 
vector  functions.  An  advantage  of  this  formulation  is  that  gradients  of  scalar 
functions  (such  as  the  pressure  term  in  the  Navier-Stokes  equations)  are 
orthogonal  to  the  divergence-free  space34  and  so  a  result  of  a  Galerkin  for¬ 
mulation  is  that  pressure  is  eliminated  entirely  from  the  problem.  A  compu¬ 
tational  method  of  this  type  was  used  be  Leonard  and  Wray35  in  their  study 
of  the  stability  of  pipe  flows. 

Throughout  our  study  the  underlying  forcing  frequency  was  taken  to  be 
n=4  and  so  that  the  critical  value  of  the  bifurcation  parameter  ft  is 

Qc  =  4v/2.  (3.5) 


Furthermore,  we  shall  be  interested  in  the  long-term,  asymptotic  behavior  of 
the  flows  developing  from  the  primary  flow 


rr  s>n  4  y 
u  =  U  =  — e. 


which  is  ;t  fixed  point  solution  of  (2.8,0).  To  facilitate  this  study  wo  use  a. 
small  fixed  perturbation  of  the  primary  flow  including  all  positive  and  nega¬ 
tive  wavenumbers  as  the  initial  conditions  for  the  calculations.  In  particular. 
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we  set 


1  If)-4 

uo(x,0)  =  —  sin4 2/e*  + -7^7  Y1  ((1  +  *)Bm(x)  -  (1  -  0B-m(*)]  (3.7) 


mi  >0 
m2  >0 


If  we  write 


u0  =  (u0,u0) 


(3.8) 


then  for  example  uQ  has  the  form 

l 


Cl 


sin  Ay  +  10  4  V'  -(cos^n^x  +  m2y)  -  sin(mix  +  m2y ) 

ni 

.  mi>0  1  1 

\  mj50 


\ 

/ 

(3.9) 


The  computations  were  performed  for  the  range 


Cl/Clc  <  12.5 


(3.10) 


After  allowing  a  considerable  time  for  the  calculation  to  settle  down  (on 
average  5000  time  units),  we  examine  the  resulting  flow  on  the  attractor.  It 
was  found  out  that  a  32  X  32  spatial  grid  was  adequate  for  the  resolution 
at  the  ratios  (3.10).  Two  types  of  dealiasing  were  used  during  the  course  of 
the  simulation  as  were  grids  of  twice  the  resolution.  Neither  were  found  to 
significantly  alter  the  '•'•suits. 
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4  Tools  of  the  Analysis 

Here  we  briefly  comment  on  the  tools  used  in  the  analysis  of  the  di¬ 
rect  numerical  simulation  of  the  Kolmogorov  flow.  The  tools  needed  for  the 
investigation  can  be  grouped  according  to  their  function. 

4.1  Temporal  Power  Spectra 

We  need  some  mathematical  tools  to  distinguish  between  chaotic  and 
quasiperiodic  regimes  of  the  flow.  The  temporal  power  spectra  (PS)  are 
indispensable  for  this  purpose30. 

Let  PS(~,(t))  denote  the  temporal  PS  of  the  time  signal  It  follows 

from  equation  (3.3)  that  there  are  322  =  1024  individual  time  histories  avail¬ 
able  for  the  use  in  the  temporal  PS.  For  the  problem  at  hand,  we  define 
the  energy  PS  box  averaged  per  7/? ,  wavenumber  in  ^-direction  as  (modulo 
a  multiplicative  normalization) 

/W(u)(”-'!  =  E/,^^(^1.rn2)(0))  +  rs(hn("i,„l.m,)(i))).  (1.1) 

in  2 

In  other  words,  we  take  the  temporal  PS  of  all  the  time  histories  result- 


ing  from  t  he  pro  jo 

ct  ions  out  o  b, 

asis  function  (3.2)  which 

indud 

le  Hi,  in  the 

wavenumber  pair  in  -  (m1.in •_> ) 

and  average  over  t  hem 

1'lns  i 

S  III  keeping 

with  i  lie  definit  ion 

of  the  PS  of  ; 

i  time  st  at  i<  ilia  ry  signa  1  a- 

-  beui!’ 

i  the  |(  inner 

transform  ol  the  autocorrelation.  In  the  calculat  ions  we  frequent  Iv  sei  m,  to 
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zero  and  use  the  resulting  temporal  power  spectrum  to  distinguish  between 
chaotic  and  quasiperiodic  regimes  of  the  flow. 

For  later  reference  we  remark  now  that  since  the  PS  is  the  Fourier  trans¬ 
form  of  the  autocoi relation  it  follows  that  as  the  frequency  |  /  |— >  oo  the  PS 
vanishes  faster  than  any  inverse  power  of  |  /  |.  This  is  a  simple  consequence 
of  the  fact  that  in  any  numerical  simulation  we  are  integrating  a  system  of 
ordinary  differential  equations  with  smooth  direction  fields  and  hence  the 
solution  is  Cco  in  the  time  variable. 

4.2  Poincare  Sections 

We  use  Poincare  sections  to  study  geometrical  structure  of  the  attractor 
and  to  qualitatively  discuss  its  degree  of  complication.  Methods  for  the 
calculation  of  these  sections  have  been  given  by  Keefe5  and  others.  The 
phase  space  employed  in  the  calculation  of  the  Poincare  sections  is  spanned 
by  the  spatial  Fourier  modes  of  the  solution  u(x./) 

U(x,t)  =  XP''»(0Dm(x)  (1.2) 

111 

where  tim  and  Bm  are  complex. 

I  he  choice  of  the  sect ’Killing  plane  for  the  Kolmogorov  Mow  becomes  very 
important,  lor  it  determines  which  attractor  features  are  covered  bv  the 
st.ndv.  Here,  we  arc  mainiv  concern  with  the  visual  structure  of  the  at 


tractor  it  self.  ('PI  expend 1 1  ures  to  aecumulal < •  a  la rge  enough  collect  ion  of 


points  for  a  highly  resolved  Poincare  section  become  considerable.  Therefore 
a  need  arises  for  the  use  of  the  discrete  symmetry  group  of  the  flow  to  extend 
the  available  data.  This  in  turn  suggests  that  some  care  should  be  exercised 
in  the  choice  of  a  sectioning  plane.  In  particular,  an  acceptable  choice  is  a 
plane  which  is  invariant  under  the  symmetry  group.  It  is  also  a  good  practice 
to  choose  as  a  section,  a  plane  for  which  the  flow  projections  maximize  the 
energy.  With  that  in  mind,  we  take  the  forcing  function  plane 


Irn(m  o,4))  =  D 

(4-3) 

Im{u(o,4))  >  0 

(4-4) 

as  our  invariant  under  the  action  of  the  symmetry  group  section  plane  for  the 
Poincare  sections  used  later.  The  choice  of  the  constant  D  controls  apparent 
positioning  of  the  hyperplane  in  the  phase  space  and  is  dependent  upon  Q. 
In  actual  calculations,  we  consider  a  series  [Im{ul0A)(t ,))];1i  of  snapshots  of 
th<>  variable  /m(?/(0i., )(/)).  We  then  set.  I)  to  its  average,  i.e. 

1  N 

D  =  -Tj  (-i..r)) 

'  1  =  1 

l  he  Poincare  sections  shown  in  the  later  sections  plot  /tV(i/(()1))  vs  /m(u((U)). 

for  further  reference,  we  note  hero  representation  of  the  symmetrv  group 
(2.12.1-1)  in  terms  of  the  spectral  components  of  u(x./).  If  u(x,/)  is  rep 
resented  as  in  (1.2)  then  the  glide-reflection  r/u(x./)  (2.12)  can  be  written 

as 


12 


gu(x,t)  =  ^2-u{_mum2)(t)el^K B,„(x)  =  J2  Vu,-m2)(0e'^Bm(x) 

in  m 

(4.6) 

while  the  rotation  Ru  (2.13)  is 

Ru(x,  t)  =  ^2  u_m(0Bm(x)  =  -um(0Bm(x)  (4.7) 

in  m 

It  then  follows  from  (4.5)  and  (4.6-7)  that  our  sectioning  plane  is  invariant 
under  the  symmetry  group  of  the  flow. 

From  the  examination  of  the  equation  (4.7),  it  is  apparent  that  the  ro¬ 
tation  corresponds  to  the  reflection  around  imaginary  axis  for  any  variable 
um(t).  On  the  other  hand  glide-reflection  action  (4.6)  of  the  symmetry  group 
on  u(o,i)(0  gives 

—  H(o,i)(0r'  *  (4-8) 

This  corresponds  to  (j7r  rotation  in  /?e(?qo,i)(0)  —  /?,,(w(o,i)(0)  plane.  Thus 
the  discrete  symmetry  group  of  the  flow  induces  an  extra  symmetries  onto 
T/(o,i ) ( )  plane.  They  correspond  to  the  reflection  around  ?  =  0.  1,2,3 

axis.  These  symmetries  figure  prominently  in  the  discussion  of  the  results  of 
our  numerical  simulat  ion. 
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4.3  Lyapunov  Spectrum 


All  of  the  mathematical  tools  mentioned  above  can  indicate  the  presence 
of  chaos  in  the  dynamical  system,  hut  do  not  directly  address  the  level  of 
randomness  present,  or  how  chaotic  a  solution  is. 

We  compute  the  Lyapunov  exponents  for  the  two  dimensional  Navier- 
Stokes  equations  (2.8,9)  to  study  the  dimension  of  the  attracting  set  of  the 
chaotic  flow.  We  write  the  governing  system  in  symbolic  form  as 

jjf  =  -Flu).  (19) 

and  denote  a  reference  or  fiducial  trajectory  by  u".  If  an  initially  nearby 
trajectories  is  written  as 

11  =  11°  +  tv  (1.10) 


then  in  the  limit  (  -*  0,  v  satisfies  the  variational  equation, 


-V  =  £(u> 


(1  11) 


where 


£K) 


^(u") 

<Su° 


is  easily  obtained  through  t  he  lineari/at  ion. 


(1.12) 


1  he  largest  Lyapunov  exponent  is  defined  as  t  he  maximum  exponential 
rate  at  which  nearby  trajectories  diverge 


A  i 


(1.1  t) 


I  i 


Thus  the  idea  of  a  Lyapunov  exponent  can  be  regarded  as  a  generalization 
of  the  stability  exponent  of  an  equilibrium  solution,  and  finding  a  positive 
Lyapunov  exponent  is  an  unambiguous  signature  of  a  chaos. 

In  a  practical  sense  we  expect  all  initial  data,  v(0),  to  lead  to  Aj,  i.e., 
as  time  evolves  v  will  become  parallel  to  the  direction  of  greatest  growth. 
Thus  to  obtain  the  next  greatest  growth  rate,  and  so  forth  can  be  a  touchy 
problem.  Algorithms  for  calculating  these  has  been  given  by  Bennetin  et 
al3‘,  Shimeda  and  Nagashima38  and  Wolf  et  al39,  and  are  further  explained 
by  Fitzmaurice40.  The  essence  of  the  trick  lies  in  computing  volumes  of 
infinitesimal  ellipsoids.  In  fact  if  en(r )  denotes  the  volume  of  an  n-ellipsoid 
it  then  follows  that 


Ai  +  A2  +  . . .  +  An  =  liin  —  In  cn(r).  (‘1.14) 

r  — *oo  7- 

Thus  after  (4.13)  one  calculates  Ai  +  A2,  then  A!  +  A2  +  A3,  etc  and  thereby 
the  spectrum  A 1 ,  A2, . . .. 

Since  the  system  in  question  is  dissipative  we  know  that  the  sum  (4.14) 
will  become  negative  for  n  large  enough.  And  since  it  is  chaotic  it  must  be 
positive  for  n  =  1.  Therefore  this  sum  must  cross  zero  at  some  point  as  one 
increases  11.  One  may  then  use  the  formula  of  Kaplan  and  Yorke22 


to  estimate  the  dimension,  d /y,  of  the  attractor  where  N  is  the  largest  integer 
for  which  the  sum  is  positive.  Constantin  et  al41  have  shown  that  the  Kaplan- 
Yorke  formula  gives  an  upper  hound  for  the  Hausdorf  dimension. 


Hi 


5  Results 


In  this  section  we  present  the  results  of  our  direct  numerical  simulation 
of  the  Kolmogorov  flow.  We  have  found  that  there  are  two  main  regimes 
of  the  flow  corresponding  to  different  ranges  of  the  bifurcation  parameter 
fi.  Tentatively  these  are  termed  as  being  the  regimes  of  small  and  large 
scale  structures.  Explanation  for  this  terminology  will  be  evident  from  the 
streamline  plots  to  be  discussed. 

5.1  Small  Scale  Structures 

o  <  n/nc  <  i 

For  this  range  the  asymptotic  solution  of  the  N-S  equations  (2.8,9)  is 
the  fixed  point  (3.6).  Streamlines  and  vorticity  plots  for  the  parallel  shear 
flow  in  horizontal  direction  are  shown  in  Figure  1. 

1  <  Q/Qc  <  1.97 

As  we  increase  Q,  the  horizontal  shear  flow  becomes  unstable  and  the 
flow  pattern  changes.  At  U/Qc  =  1,  the  parallel  shear  flow  bifurcates  to 
a  new  solution.  In  this  range  the  motion  is  steady  and  purely  cellular  and 
illust  rations  of  t  his  are  shown  in  Figures  2  and  3.  It  is  interesting  to  note  that 
at  the  lower  range  of  the  parameter  U  the  flow  is  almost  horizontal  (Figure 
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2a)  while  at  higher  values  of  Q  in  this  range  the  motion  is  mainly  in  the 
vertical  direction  (Figure  3a)  with  well  defined  cells.  In  the  first  instance  the 
saddle  points  of  the  separatrices  are  formed  horizontally  while  in  the  latter 
case  they  are  joined  vertically.  This  is  in  an  agreement  with  the  experiments 
of  Bondarenko,  et  al42  for  a  magnetohydrodynamic  problem  under  similar 
forcing. 

1.97  <  f yne  <  2.22 

In  this  range  the  flow  still  have  the  steady  cellular  structure  but  the 
number  of  cells  in  the  horizontal  direction  have  doubled.  The  new  flow 
patterns  are  depicted  in  Figure  4.  Explanation  of  this  phenomena  relies  on 
the  fact  that  the  horizontal  wavenumber  of  the  fastest  growing  mode  changes 
from  1  to  2,  also  reported  by  Green16.  A  very  important  feature  of  the  flow 
in  this  and  previous  ranges  is  that  it  is  invariant  under  the  discrete  symmetry 
group  transformations.  By  invariance  we  mean  that  for  any  transformation 
(!  in  the  discrete  symmetry  group  of  the  flow  there  is  0  <  /  <  2tt  such  that 

1)0 vi  ~  u  (5.1) 

where  7)  is  a  translation  by  l  in  .r-direetion  and  u  is  a  solution  of  N-S  equa¬ 
tions  (2.X. 9).  In  words,  we  say  that  the  flow  is  invariant  under  the  discrete 
symmetry  group  if.  after  application  of  any  trauslormat  ion  from  the  discrete 
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symmetry  group  to  the  flow,  we  can  translate  the  transformed  flow  in  liori 
zontal  direction  in  such  a  way  that  it  coincides  with  the  original  flow. 


2.2  <  n/Uc  <  2.53 

As  f?  increases,  small  amplitude  perturbations  at  the  saddle  points 
will  lead  to  the  destabilization  and  eventual  breakdown  of  the  steady  cellular 
pattern.  At  Q/flc  —  2.2  another  bifurcation  in  parameter  space  takes  place. 
For  me  range  under  discussion  there  is  a  gradual  breakdown  in  the  symmetry 
of  the  flow.  The  new  flow  patterns  are  depicted  in  Figure  5.  It  is  perhaps 
remarkable  that  these  cellular  structures  are  still  steady. 

3.53  <  n/Qc  <  3.6 

At  Q/S}c  —  3.53  the  statistically  stationary  flow  is  no  longer  steady. 
Now,  the  flow  spends  most  of  the  time  in  the  steady  cellular  state  shown 
in  Figure  5.  But  from  time  to  time,  the  flow  is  intermittent  ly  disturbed  by 
a  random  horizontal  shifts  orruring  on  a  very  rapid  time  seal*'.  The  time 
between  bursts  is  very  long  and  is  a  decreasing  function  as  the  parameter 
S2  reaches  the  upper  end  of  the  range.  Figure  6  shows  the  time  series  <>f 
Ilc(v^)  i j ).  Here,  the  plateaus  correspond  to  a  steady  cellular  stale.  After  the 
horizontal  shift,  the  cells  oscillate  with  a  small  amplitude  for  a  short  while. 
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and  then  settle  down  very  rapidly  to  a  steady  pattern  before  undergoing  the 
next  horizontal  shift.  This  state  is  chaotic  but  only  over  a  very  long  time 
scale,  i.e.  the  Lyapunov  spectrum  has  some  very  small  positive  exponent. 

3.6  <  «/ftc  <  4.47 

In  this  range  the  nature  of  the  flow  changes  once  again.  We  now  find  a 
periodic  oscillations  occuring  in  the  cellular  structures  as  depicted  in  Figure 
5.  Also,  the  solution  develops  a  slow  drift  in  horizontal  direction,  but  without 
a  net  transport  off  mass.  The  drift  may  be  regarded  as  the  wave  motion.  It. 
also  follows  from  symmetry  that  the  drift  wave  can  move  in  either  direction. 
We  conclude  from  Figure  7  which  shows  the  temporal  power  spectrum  that 
there  is  a  principal  frequency,  the  large  peak,  and  subharmonic  frequencies 
of  it  present  in  the  flow.  These  correspond  to  the  oscillations  of  the  pattern 
around  larger  vortices  and  inner  oscillations  of  the  vortices  themselves.  It  is 
interesting  to  note  that  we  have  found  no  regime  of  the  flow  when  only  the 
principal  frequency  is  present. 

4.47  <  n/Qc  <  4.58 


Now,  the  flow  enters  its  second  chaotic  regime.  Figure  8  depicts  a 
typical  time  series  of  the  signal  lfc( ?/ (0  i ) )  in  this  range,  from  which  we  con 
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elude  that  the  chaos  is  intermittent.  Examining  the  power  spectrum  shown 
in  Figure  9,  we  see  that  the  principal  and  subharmonic  frequencies  are  promi¬ 
nently  present  in  what  is  a  broad  hand  spectrum.  We  did  not  attempt  to 
calculate  the  Kaplan- Yorkc  dimension  for  the  dynamical  system  in  this  range 
of  fl  because  of  the  very  long  time  interval  between  the  bursts  of  activity  in 
the  flow.  A  most  remarkable  feature  of  the  flow  in  this  range  is  that  the 
chaotic  attractor  is  not  ergodicl  The  dynamical  system  has  several  chaotic 
attractors  each  with  its  own  basin  of  attraction  lying  quite  close  to  the  fixed 
point  (3.6).  On  the  basis  of  extremely  long  calculations  we  conclude  that 
these  attractors  are  not  connected.  Figure  10  shows  the  Poincare  section 
depicting  cuts  of  the  chaotic  attractors.  Two  darker  points  in  the  blown  up 
section  Figure  10b  correspond  to  the  punctures  produced  by  the  solution  in 
the  laminar  phase  of  the  intermit  ten  cy.  In  actual  calculations  we  obtain  only 
one  of  the  attractors  (depending  on  the  chosen  initial  data)  ,  and  the  rest 
are  obtained  by  using  the  discrete  symmetry  group.  Alternatively,  we  can 
apply  the  symmetry  group  to  the  initial  data  and  thus  generate  each  of  the 
attractors.  Another  peculiarity  of  the  flow  is  that,  there  are  only  eight  sym¬ 
metric  attractors  instead  of  sixteen  as  predicted  by  the  discrete  symmetry 
group.  I'h is  happens  because  the  action  of  five  glide-reflections  followed  by 
the  rotation  by  7T  applied  to  the  solution  maps  the  attractor  onto  itself.  This 
induces  an  extra  symmetry  into  the  attractor  corresponding  to  the  reflection 


around  the  axis  denoted  by  the  dashed  line  in  Figure  10b,  as  it  is  discussed 
in  the  earlier  part  of  this  paper  dealing  with  the  Poincare  sections. 

4.58  <  fi/nr  <  4.65 

As  Cl  increase,  different  symmetric  attractors  start  intersecting  each 
other  and  now  the  flow  is  intermittenly  chaotic  with  one  ergodic  attractor 
consisting  of  eight  connected  symmetric  components.  Figure  11  shows  the 
time  series  of  /<,e(tq0,i))-  The  evidence  of  the  intermittent  bursts  is  clearly 
present  and  sudden  jumps  correspond  to  the  solution  jumping  to  a  different 
component  of  the  attractor.  We  present  the  temporal  power  spectrum  in 
Figure  12.  This  new  state  is  more  chaotic  than  the  previous  case,  but  still 
there  is  evidence  of  the  principal  frequency  from  the  lower  parameter  values. 
'Die  Poincare  section  is  shown  in  Figure  13  (“butterflies”  section).  Darker 
regions  in  the  blown  up  section  of  the  attractor  in  Figure  13b  correspond  to 
the  points  where  the  solution  spends  most  of  the  time  between  intermittent 
bursts.  Moreover,  comparing  Figures  13  and  10.  we  notice  that  the  main 
features  of  t  he  non-ergodie  attractor  are  still  present. 


1  65  •  /  <)/(>  <  4.71 


Another  remarkable  feature  of  the  Kolmogorov'  How  is.  that  counter 
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to  intuition,  the  parameter  space  fl  contains  windows  of  relaminarization 
between  different  chaotic  ranges  of  the  flow.  This  range  is  an  example  of 
one  of  such  windows.  Here,  the  flow  returns  to  a  periodic  regime  with  a 
principal  and  subharmonic  frequencies  in  the  power  spectrum.  Again,  ihero 
is  no  regime  with  just  one  frequency.  This  flow  is  similar  to  the  range  3.6  < 
fl/flc  <  4.47.  It  also,  develops  a  slow  steady  drift  in  the  horizontal  direction. 
We  also  note  that  this  drift  appears  to  vanish  in  the  chaotic  windows  of  the 
flow. 

At  fl  =  4.7'2flc  we  see  the  formation  of  large  scale  structures. 

5.2  Large  Scale  Structures 

Although  we  are  mainly  interested  in  the  asymptotic  behavior  of  the 
Kolmogorov  flow,  it  is  of  interest  to  examine  how  the  large  scale  structures 
are  formed.  The  process  is  basically  the  same  for  all  values  fl  >  4.72Q,.. 

We  start  with  a  small  perturbation  of  the  fixed  point  flow  (3.6).  This  state 
persists  for  roughly  the  time  scale  predicted  by  the  linear  stability  analysis. 
Then  on  a  very  rapid  time  scale  the  flow  changes  to  the  unsteady  cellular 
pattern  of  the  sort  illustrated  in  Figure  5.  T  his  state  with  many  small  cells 
persists  for  a  very  short  time  and  then  the  flow  plunges  into  chaos.  The  new 
state  is  vigorously  chaotic,  characterized  bv  t  'ie  rapid  movements  of  a  small 
vortical  structures,  a  snapshot  of  which  is  shown  in  Figure  1  1.  |?v  a  rough 


count  it  is  seen  that  the  same  number  of  vortical  structures  appears  in  Figure 
5  and  14.  We  refer  to  this  state  as  metastable  or  transient  chaos  following  the 
terminology  of  Grebogi,  Ott  and  Yorke43.  Of  all  the  chaotic  states  observed 
in  our  numerical  experiments,  this  is  the  only  state  which  can  be  called 
chaotic  both  spatially  and  temporally  and  hence  is  turbulent  according  to 
common  convention.  At  the  lower  range  of  fl  this  metastable  chaos  holds  for 
a  relatively  long  time,  while  at  a  higher  values  of  the  bifurcation  parameter 
it  can  be  extremely  short. 

In  any  case  the  phase  point  finds  an  escape  hole  in  the  metastable  attrac¬ 
tor  and  we  eventually  sec  the  formation  of  a  large  scale  structures  such  as 
those  depicted  in  the  Figure  15.  The  new  attractor  is  stable.  Considerable 
vortex  merging  has  taken  place  in  passing  from  Figure  14  to  Figure  15.  The 
two  large  structures  shown  in  the  streamline  plot  of  Figure  15a  are  some¬ 
what  misleading.  On  a  closer  examination  of  the  vorticity  plot,  wre  see  that 
we  have  roughly  the  same  number  of  distinct  vortices  as  before  and  that  the 
large  structures  in  the  streamline  plot  consist  of  a  number  of  a  counterrotat¬ 
ing  vortices.  These  new  large  scale  structures  are  completely  different  from 
the  small  scale  structures  which  we  dealt  in  the  previous  subsection,  exam¬ 
ining  the  vorticity  plot,  shown  in  f  igure  5b.  we  notice  that  at  the  small  scales 
the  vortices  arranged  in  a  lattice.  Thus,  the  unsteady  small  scale  st  nut  ures 
flow  can  be  considered  as  oscillations  or  '  breathing  of  the  lattice  structure. 


In  the  metastable  regime  this  lattice  structure  is  broken,  and  we  have  a  soup 
of  vortices  moving  randomly  around.  Later  the  vortices  merge  together  and 
arrange  themselves  on  a  curved  line  with  two  vortices  containing  most  of  the 
circulation.  These  vortices  correspond  to  the  eves  of  the  large  scale  struc¬ 
tures  depicted  in  the  streamline  plot  of  Figure  15a.  It  is  of  interest  to  note 
here  that  the  vertical  position  of  the  eyes  of  the  large  scale  structures  are  at 
hyperbolic  points  of  the  forcing  functions  when  the  flow  is  not  chaotic. 

We  will  now  describe  in  detail  the  bifurcations  occuring  in  the  large  scale 
structure  regime  of  the  flow. 

4.72  <  n/nc  <  5.9 

In  this  range  the  flow  is  periodic  with  one  frequency.  In  addition,  there 
is  a  horizontal  drift  present  in  the  solution  as  also  noted  by  She26.  As 
approaches  5.9  both  the  frequency  of  oscillations  and  drift  speed  approach 
zero.  Figure  16  shows  the  temporal  power  spectrum. 

Q  =  5.9QC 

At  il  =  the  flow  reaches  the  steadv.  spatiallv  disorganized  state 

depicted  in  Figure  15.  The  flow  Is  chaotic  only  in  the  spatial  domain! 
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5.9  <  n/nc  <  9.47 


As  we  start  increasing  Q.  further,  the  flow  develops  both  periodic  oscil 
lations  and  a  slow  drift  in  the  horizontal  direction.  The  direction  of  the  drift 
changes  according  to  the  particular  attractor  the  flow  settles  in.  Moreover, 
as  Q  increases,  so  does  the  drift  speed  and  frequency  of  the  oscillations.  The 
periodic  oscillations  can  be  described  as  the  “breathing”  in  the  inner  shells 
of  the  large  scale  structures.  At  fi  =  9.042QC  a  period  doubling  bifurcation 
takes  place  in  the  oscillations  of  the  inner  structures.  Now  the  principal  and 
subharmonic  frequencies  correspond  to  the  breathing  of  the  large  scale  struc¬ 
ture  itself  and  the  oscillations  of  the  inner  shells.  This  regime  persists  up 
to  Q  =  9.450c,  where  another  period  doubling  bifurcation  in  the  oscillations 
takes  place.  We  can  see  clearly  the  principal  and  subharmonic  frequences  in 
Figure  17  showing  the  temporal  power  spectrum  of  the  (low. 

9.47  <  ^  9.f»fi 

At  ff  =  9.47U,-  the  flow  <-nters  another  chaotic  regime  which  continues 
until  V.  reaches  9.(iGfb  .  We  again  find  the  remarkable  feature  in  this  win¬ 
dow  that  there  is  chaos  but  that  the  attractois  are  non-t  rqnthr]  Here,  the 
chaos  is  present  in  t  he  irregular  breat  hing  like  u-,<  i 1 1  a t  ions  of  t  he  inner  U  rue 
I  ii res  which  induces  random  < >sr dial  ions  of  la rge  s<  ale  si  met  u res.  We  observe 
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in  Figure  18,  which  shows  the  temporal  power  spectrum,  that  the  driving 
frequency  and  its  subharmonics  still  dominate  the  motion.  Horizontal  drift 
diminishes  to  a  trickle  and  sometimes  changes  directions.  On  average,  there 
is  no  drift  at  all.  Figure  19  depicts  the  Poincare  section  of  the  attractor 
showing  eight  symmetric  disconnected  attractors.  The  slopes  of  the  lines  in 
Figure  20  correspond  to  the  values  of  Lyapunov  exponents  and  the  calculated 
Kaplan-Yorke  dimension  is  7.9. 

9.66  <  Q/Qc  <  10.805 

At  ft  =  9.66HC  the  flow  enters  another  laminar  window  in  the  pa¬ 
rameter  space  (l.  Here,  the  solution  slowly  drifts  in  the  horizontal  direction 
and  has  the  principal  and  subharmonic  frequencies.  At  H/Hc  =  10.4  a  pe¬ 
riod  doubling  bifurcation  takes  place.  Figure  21  shows  the  temporal  power 
spectrum  at  H  =  10.7HC. 

10.805  <  n/nc  <  12.5 

As  we  increase  H  even  further,  the  flow  enters  another  chaotic  regime  at 
H  =  10.805H,'-  The  flow  now  appears  to  bo  governed  by  an  orgodio  attractor 
which  consists  of  eight  symmetrical  components  and  (In*  solution  intermit 
tently  jumps  from  one  component  to  another.  There  is  a  very  long  time 


delay  between  successive  jumps  at  0  slightly  greater  than  the  10.805,  but 
it  rapidly  decreases  as  0  increases.  The  driving  frequency  and  its  subhar¬ 
monics  dominate  the  power  spectrum  shown  in  Figure  22.  An  examination 
of  Figure  23  showing  the  Poincare  section  of  the  chaotic  attractor,  indicates 
that  the  components  of  the  attractor  have  a  fractal  structure.  The  calculated 
Kaplan- Yorke  dimension  is  8.9.  Another  interesting  point  is  that  the  direc¬ 
tion  of  the  drift  changes  with  each  successive  jump  to  a  different  component 
of  the  attractor. 

As  we  continue  to  increase  Cl,  the  different  components  of  the  chaotic 
attractor  merge  together  and  become  indistinguishable.  Nevertheless,  the 
driving  and  subharmonic  frequencies  are  still  dominate  the  temporal  power 
spectrum  shown  in  Figure  24  for  0  =  12.50,-.  Now,  the  Kaplan-Yorke  di¬ 
mension  increases  to  10.8.  Thus  we  might  hope  to  describe  the  system  in 
terms  of  a  dynamical  system  having  many  fewer  than  roughly  1000  degrees 
of  freedom  that  we  have  been  allowing  in  our  spectral  simulations. 

As  in  the  case  of  the  small  scale  structures,  the  average  drift  speed  di¬ 
minishes  to  zero  for  the  chaotic  regimes  and  the  drift  is  not  present  at  all  for 
0  >  12.50,.. 

I  his  conclude:;  our  detailed  description  of  t lie  behavior  of  the  fixed  point 
(3.0)  of  the  Kolmogorov  flow  as  a  function  of  the  bifurcation  parameter  O. 


28 


6  Summary 


In  this  paper  vve  investigated  in  some  detail  a  simple  fluid  flow  governed  by 
the  incompressible  two-dimensional  Navicr-Stokes  equations  with  a  spatially 
periodic  time  independent  forcing.  Two  main  regimes  of  the  flow  have  been 
observed:  small  and  large  scale  structure  regimes  corresponding  to  different 
ranges  of  Reynolds  number  Re.  A  sequence  of  bifurcations  takes  place  in 
each  of  the  regimes  of  the  flow. 

The  first  three  bifurcations  taking  place  in  the  small  scale  structure  regime 
are  steady:  the  first  corresponds  to  the  appearance  of  the  cellular  structures, 
the  second  doubles  the  horizontal  wavenumber  of  the  cellular  structures, 
while  the  third  breaks  the  symmetry  of  the  flow.  In  addition,  three  chaotic 
windows  are  discovered  in  the  small  scale  structures  regime.  The  first  window 
is  chaotic  over  a  very  long  time  scale.  The  second  window  is  an  intermittency 
with  a  remarkable  feature  that  the  chaotic  flow  has  non-ergodic  attractors. 
Finally,  the  last  chaotic  window  shows  intermittency  with  one  connected 
chaotic  attractor  consisting  of  eight  symmetric  components  (“butterflies”). 
There  arc  laminar  windows  between  the  first  and  the  second  chaotic  regimes, 
and  after  the  last  chaotic  regime. 

The  large  scale  structure  regime  contains  a  number  of  quasi -periodic  and 
chaotic  windows.  Before  the  formation  of  the  large  scale  structures,  the  flow 


spends  considerable  time  in  metastable  chaos  which  approaches  true  turbulent 
flow  for  the  range  of  Re  under  the  discussion.  A  spatially  disordered,  but 
temporally  steady  state  was  observed  between  first  two  laminar  windows. 
Two  period-doublings  take  place  in  the  second  laminar  window.  The  chaotic 
regime  which  follows  has  the  same  remarkable  property  of  chaotic  but  non- 
ergodic  attractors  as  is  the  case  with  the  small  scale  structure  regime.  In 
addition,  there  is  a  window  of  relaminarization  with  a  quasi-periodic  flow 
between  the  two  chaotic  regimes. 

In  the  future,  we  plan  to  use  the  proper-orthogonal  decomposition44  to 
extract  the  coherent  structures  of  the  flow  and  study  the  resulting  low  order 
dynamical  systems45. 
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1.  Parallel  shear  layer  flow  corresponding  to  the  fixed  point  solution  at 
=  0.5fic.  a)  streamlines. 
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2.  Steady  cellular  flow  found  at  Q  =  1.12fie  with  the  motion  mainly  in 
horizontal  direction,  a)  streamlines. 
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2.  Steady  cellular  flow  found  at  ft  = 
horizontal  direction,  b)  vorticity. 


1.12ne  with  the  motion  mainly  in 


41 


3.  Steady  cellular  flow  found  at  fl  =  1  -5f2c  with  the  motion  mainly  in  ver¬ 
tical  direction  and  the  fastest  growing  perturbation  with  wavenumber 
1  in  horizontal  direction,  a)  streamlines. 
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3.  Steady  cellular  flow  found  at  fi  =  1.5nc  with  the  motion  mainly  in  ver¬ 
tical  direction  and  the  fastest  growing  perturbation  with  wavenumber 
1  in  horizontal  direction.  b)vorticity. 


tical  direction  and  the  fastest  growing  perturbation  with  wavenumber 
2  in  horizontal  direction,  a)  streamlines. 


4.  Steady  cellular  flow  found  at  fl  =  2fle  with  the  motion  mainly  in  ver¬ 
tical  direction  and  the  fastest  growing  perturbation  with  wavenumber 
2  in  horizontal  direction,  b)  vorticity. 
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6.  Time  series  of  J2e(u(o,i))(<)  at  I)  —  3.56fle> 


steady  cellular  flow. 
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3.0E4 


4.0E4 


Plateaus  correspond  to  a 


10.  a)  Poincare  section  of  the  chaotic  attractors  at  fl  =  4.48ftc.  Attractors 


are  not  ergodic. 


SB 


13.  a)  Poincare  section  of  the  chaotic  attractor  at  ft  =  4.63QC.  Flow  inter¬ 
mittently  jumps  froih  one  component  of  the  attractor  to  another. 


06 


13.  b)  One  of  the  components  of  the  attractor. 
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14.  Instantaneous  flow  pattern  found  in  metastable  chaotic  regime  at  ft 
10ftc.  b)  vorticity.  The  flow  is  both  spatially  and  temporally  chaotic, 

and  thus  may  be  called  turbulent. 
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15.  Steady  spatially  chaotic  flow 

this  is  an  example  of  large  scale  structures 


found  at  0  =  5.9SV  »)  streamlines. 


Also, 
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15.  Steady  spatially  chaotic  flow  found  at  fl  =  5.9f lc. 
this  is  an  example  of  large  scale  structures. 
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16.  Power  spectrum  of  the  velocity 
present  in  the  problem. 
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19.  a)  Poincare  sections  of  the  chaotic  attractor  at  Q  =  9.47fic.  Again, 
attractors  are  not  ergodic. 
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Frequency  /  Hz 


22.  Power  apectrum  of  the  velocity  u  at  fl  =  10.805SJ,.  The  flow  i>  chaotic, 
but  the  driving  frequency  end  its  subharmonica  are  dominant  m  the 

spectrum. 


23.  c)  Blow  ups  of  one  of  the  components  of  the  attractor.  The  fractal 


structure  of  the  attractor  is  evident  from  blow  up. 
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